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Abstract 

We present a numerical study of mixed boson-fermion systems at zero tem- 
perature in isotropic and anisotropic harmonic traps. We investigate the 
phenomenon of component separation as function of the strength of the inter- 
particle interaction. While solving a Gross-Pitaevskii mean field equation for 
the boson distribution in the trap, we utilize two different methods to extract 
the density profile of the fermion component; a semiclassical Thomas- Fermi 
approximation and a quantum mechanical Slater determinant Schrodinger 
equation. 
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I. INTRODUCTION 



Since the recent experimental realization of Bose-Einstein condesation in dilute gases of 
rubidium sodium @,||, lithium |7|], and hydrogen a great deal of interest in Bose 

condensed systems has concentrated on the topic of multi-component condensates. This 
field was stimulated by the succesful demonstration of overlapping condensates in different 
spin states of rubidium in a magnetic trap an d of sodium in an optical trap [ I I 



the (binary) mixtures being produced either by sympathetic cooling, which involves one 
species being cooled to below the transition temperature only through thermal contact 
with an already condensed Bose gas, or by radiative transitions out of a single component 
condensate. Since then a host of experiments has been conducted on systems with two 
condensates, exploring both the dynamics of component separation |]12;] , and measuring the 
relative quantum phase of the two Bose-Einstein condensates [[L3|]. Most of the theoretical 
work concerning multi-component condensates [H-^] has been devoted to systems of two 



Bose condensates. However, other systems are of fundamental interest, one of these being 
a Bose condensate with fermionic impurities, a system reminiscent of superfluid 3 He- 4 He 
mixtures. In particular the possibilty of sympathetic cooling of fermionic isotopes has been 
predicted in both 6 Li- 7 Li 39 K- 40 K, and 41 K- 40 K [p5|1 . Magneto-optical trapping of the 



fermionic potassium isotope 40 K has been reported 

The boson-fermion mixture was discussed in a previous paper [^] within the Thomas- 



Fermi approximation, which amounts to neglecting the kinetic energy of the bosons, and to 
apply a semi-classical filling of phase space of the fermions. For the bosons, this is a valid 



approximation in the limit of strong interactions or large particle numbers, see |28| . In this 
paper we present a numerical analysis of the system, incorporating the correct operator form 
of the kinetic energy of the particles. 

The paper is structured as follows. In Sec. II we study in detail the case of an isotropic 
external potential and we develop both the Thomas- Fermi approximation and the full quan- 
tum mechanical description of the fermions. The numerical procedure is briefly introduced. 
In Sec. Ill the case of the anisotropic harmonic oscillator trap is outlined within the Thomas- 
Fermi approximation for the fermions. In Sec. IV we present our quantitative results for the 
isotropic and anisotropic trapping potentials, demonstrating the accuracy of the predictions 



made in [27|, and addressing the issue of symmetry breaking in elongated traps. Sec. V 
summarizes the main results. 

Throughout, we assume that the bosons and fermions have the same mass, M, and 
that the atoms are all trapped in the same external harmonic oscillator potential. This 
choise is of course only a convenience; all our calculations are readily generalized to differing 
experimental parameters. 

II. ISOTROPIC TRAPS 

A. Gross-Pitaevskii equation for the bosons 

In the mean field description the behavior of the single particle wavefunction ip(r), as- 
sumed to describe all Nb bosons in the gas, is governed by the Gross-Pitaevskii equa- 
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tion p9|j30|| . In the presence of fermions, this equation is modified by the addition of an 
interaction term proportional to the fermion density, uf{t) 



2M 



V 2 + V ext (r) + gN B \^(f)\ 2 + hn F (f) 



(1) 



where V ext (r) = \M{oj 2 x 2 + u 2 y 2 + uj 2 z 2 ) is the external confining potential, and \i is the 
boson chemical potential (energy per particle). The value of [i is fixed by the normalisation 
condition, / d 3 rnB{r) = Nb on the boson density Ub{t) = N b\4 ! {^)\ 2 '■ 

The low kinetic energies of the atoms permit the replacement of their short range in- 
teraction potential by a delta function potential of strength g or h. This is known as the 



pseudopotential method pi |. There is no fermion- fermion interaction in this description, 
see below. In (HD g and h thus represent the boson-boson and the boson-fermion interaction 
strengths proportional to the respective s-wave scattering lengths Pf |. 

In isotropic traps we have V ext (r) = ^Muj 2 r 2 , r being the distance from the trap center. 
By the substitution x — r " l P in (H) we obtain the radial equation 



2M dr 2 



+ 



V ext {r) + gN B 



X{r) 



+ hnpir) 



X(r) = ^x{r)- 



(2) 



In order to simplify the formalism, we rescale (§) in terms of harmonic oscillator units, that 
is 



r = dor, 

jj, = huojl, 
i[)(f) = a 3/2 ^(r), 
X(r) = % 2 x(r), 



where ao 



U/Muj is the width of the oscillator ground state. Defining 

gM ~ hM 



g 



h 



a h 2 ' aoh 2 ' 
we arrive at the simplified equation for the radial function 

_ I fx 
2 df 2 



(3) 



(4) 



+ 


\r 2 + gN B 


X{r) 


2 

+ hn F {f) 






r 





(5) 



In the remaining parts of the paper we shall omit the tilde from this equation. 

To solve the Gross-Pitaevskii equation for the bosons, we must find np(r). To this end 
we invoke two methods: A semi-classical (Thomas-Fermi) approximation and a quantum 
mechanical treatment. 
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B. Thomas- Fermi approximation for the fermions 



In the semi-classical (Thomas-Fermi) approximation the particles are assigned classical 
positions and momenta, but the effects of quantum statistics are taken into account. That is: 
The density in the occupied part of phase space is simply (2tc)~ 3 , and sums over states can 
be replaced by the corresponding integrals over for k. The fermions experience a potential 
V{f) = V ex t(r) + hnB{r) and for particle motion in such a potential it is posible to define a 
local Fermi vector kp(r) by 

so that the volume of the local Fermi sea in k space is simply 

^Kk F {rf = (27r)V(f). (7) 

In the low temperature limit, where p-wave (and higher multipole) scattering can be ne- 
glected, the supression of the s-wave scattering amplitude due to the antisymmetry of the 
many-body wavefunction implies that the spin polarized fermions constitute a noninteract- 
ing gas (for the case of an interacting Fermi gas, see ||32|| ). Hence the density of the fermionic 
component is given by 



n F (f) = {^r [Ef - V ext {r) - hn B (f)]j /(6tt 2 ). (8) 

As in the case of the bosons, where the chemical potential must be adjusted for the integral 
of the density over space to give the correct number of particles, the Fermi energy determines 
the proper normalisation; / d 3 rnp(f) = Np. For a thorough discussion of trapped fermions 
(also at T > 0), and comments on the range of validity of the Thomas- Fermi approximation 

sec 



C. Slater determinant description 

The many-body wavefunction, ty(fx . . . f^ F ), may be represented by a Slater determinant 

I n f 

^? 1 ...r Np ) = -=A]l<Pi(iV, (9) 
V i » • j=i 

where A is the antisymmetrization operator. This Slater determinant solves a stationary 
Schrodinger equation 

HV(r) = EV(r), (10) 
with a Hamiltonian that is the sum of Np independent single-particle operators 

N F 

H = J2H^ (11) 

i=i 
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The orbitals (pi(fi) solve the eigenvalue equation 

Htfiifi) = Ei<fi(fi). 

We make the substitution ip{r) = Unl ^ Y em (9, 0), where Yi m (6,< 
harmonics, and we thus obtain a radial equation for the functions u n i in harmonic oscillator 
units 

1 d 2 u ni 



(13) 

are the usual spherical 



2 dr 2 



+ 



-r" + 



£(£ + 1) 

2r 2 



+ hn B (r) 



u n e(r 



E nf u n Ar). 



(14) 



It is important to keep in mind that the radial functions must satisfy the boundary condition 
^w(O) = 0, to ensure a finite particle density at the center of the trap. 

Equation (|i"4l) can be solved once for every lvalue, thus producing the energy spectrum. 
The centrifugal term in the radial equation, implies that the fermions can be considered to 
move in an isotropic effective potential, V e ff (r) = V ext (r) + hnsir) + ■ The energy 
levels E n £ are 2£ + 1 times degenerate, and we have the fermion density given by 



n F [r\ 



E 

occupied 

states 



u n e{r) 



Yi 



r 



£(2^+1) 



ni 



\u n e(r)\' 



(15) 



since Y,m=-i \ Yem(@, 0)| 2 — (2^ + l)/47r. Once found the eigenstates are sorted after energy 
and the energy levels are filled from below with Np particles. The Fermi energy is the energy 
of the highest occupied orbital. The maximum value of the angular momentum, may be 
estimated from the Thomas- Fermi expressions (0,§) for h = by maximizing r/ci?(r), the 
maximal length of i at the point r. This yields the simple result, £ ma x ~ Ep, where the 
Fermi energy in the noninteracting limit is Ep = (6NpY^ 3 in harmonic oscillator units. To 
test our numerical calculations for fermions not interacting with the bosons (h = 0), we 
have compared our spatial density distributions with those of Schneider and Wallis []34| and 
found excellent agreement. 



D. Numerical Procedure 

We note that the solution of both Eqs. ([5]) and ( fH] ) require prior knowledge of the boson 
density %(r). To obtain the density profiles of the two components, we insert iteratively 
the density of one component into the equation for the other until a desired convergence is 
reached. 

To solve (H) for the boson density, we use the method of steepest descend, that is we 
propagate a trial function (which can be chosen initially almost arbitrarily) in imaginary 
time r, replacing nx{ r ) by — (d/dr)x{r, r). In the long time limit the propagation "filters" 
the trial function to the condensate ground state Alternative methods for solving numerically 
the Gross-Pitaevskii equation are presented in [|T9| , ^5H 37|] . 

The evaluation of the fermion density profile is done by two methods, as described above. 
In the case of the Thomas- Fermi approximation, np{r) is found by direct insertion of ub{t) 
into (|8p, searching numerically for the energy Ep giving the right number of particles. Within 
the Slater determinant method one obtains the density profile directly from (|15D , once the 



diagonalization of fll4|) has been done 
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III. ANISOTROPIC TRAPS 



In this section we treat the case of an anisotropic trapping potential with a cylindrical 
geometry (u x = ui y —u±^ u z ) as this corresponds to current experimental setups We thus 
have 

V ext = l -Mu\r 2 + l -Mu 2 z z\ (16) 

where r = \/x 2 + y 2 and z are the radial and axial coordinate respectively. We define the 
asymmetry parameter A = u z /uj±. 

As in the case of the isotropic potential we have for the bosons a non-linear Schrodinger 
equation corresponding to ([[]). By the substitution x — sft'ty we obtain the equation 



VX(r, z) = - - 



d 2 x d 2 x 
Or 2 dz 2 

|2 



^+ 1 -(r 2 + X 2 z 2 ) X (r,z) 
8r z 2 



+ 9 N B l -^)Lx(r, z) + hn F (f) X (r, z), (17) 

in harmonic oscillator units. Again the radial function has to vanish on the symmetry axis 
to remove potential problems of divergences near the origin. This boundary condition is 
implemented in our numerical procedure by imposing on the radial function a y/r dependence 
for small values of r, fitting to the value of x a ^ larger distances from the axis. As in the 
case of the isotropic trapping potential a Split-Step-Fourier technique is used to propagate 
the boson wavefunction to the ground state. An alternative method for solving the Gross- 
Pitaevskii equation in a cylindrical configuration, applying an alternating-direction implicit 
method to compute the derivatives is discussed in f38fl . 

We shall limit ourselves to the Thomas-Fermi approximation for the fermions, both out 
of necessity and convenience. Already in the spherically symmetric, effectively 1-dimensional 
case, the full quantum mechanical analysis is very time consuming and as we shall demon- 
strate in the next section, the Thomas-Fermi approximation offers the same qualitative 
features as the exact description. The fermion density is thus evaluated using equation (H) 
with the external potential ([16]) and with the boson density obtained from (|lj 



IV. RESULTS 

The main conclusion of [[^J is the prediction of a component separation under variation 
of the strength of the boson-boson and boson-fermion interaction. In the Thomas-Fermi 
approximation for both components, the density distributions solve the coupled equations 



2M 



Vext(r) + g ■ n B (r) + h ■ n F (r) = /x 
(67rV(r)) 2/3 + V ext {r) + h ■ n B {r) = E F . (18) 



In the case of Np/N B 1 the fermions may be neglected in the equation for the bosons. 
For the fermions we then obtain the simple equation 
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— (6vrV(r)) 2/3 + (1 - -)K*t(r) + -/i = (19) 

where the terms proportional with h are absent in regions with vanishing ub{t)- We may 
distinguish between 3 different types of solutions; if h < g the potential minimum of the 
fermions is located at the center of the trap, and if their number is small enough, they 
will constitute a 'core' entirely enclosed within the Bose condensate. The two quantum 
gases are truly interpenetrating. If h = g the fermions have a constant density throughout 
the Bose condensate, falling towards zero outside. If h > g the effective potential for the 
fermions is that of an inverted harmonic oscillator having a minimum at the edge of the 
Bose condensate, where the fermions localize as a 'shell' wrapped around the condensate. 



A. Isotropic trap, quantum treatment 

When we replace the Thomas-Fermi approximation by an exact description including the 
kinetic energy operator for the bosons and treating the fermions quantum mechanically, we 
expect to observe the same overall behaviour, but with minor corrections. The boson kinetic 
energy is expected to cause penetration into the fermionic component and a rounding off 
of the atomic distributions at the boundaries. Fig. [I] shows the spatial distribution of 1000 
fermions in a condensate of 10 6 bosons for different values of the boson-fermion interactions 
strength, h. The strength of the boson-boson interaction, g, is chosen to give maximal over- 
lap between the two atomic clouds. In order to have clouds of comparable size, we equate the 
Thomas-Fermi expressions for the radius of the Bose condensate (l^Nsg /AnMu 2 ) 1 ' 5 , and 
the radius of the zero temperature Fermi gas (A&Np) 1 ^ Jfi/ Mu. This gives the condition: 



g/(hua 3 ) ~ 21.1N F B '«/N B , (20) 

which for the parameters of Fig. [I] requires g = 0.015. The coupling g differs for differ- 
ent atomic species and this value is in approximate agreement with the coupling stregth 
in the MIT Na setup and we recall the possibility to achieve couplings of arbitrary 
strength by the recently demonstrated modification of the atomic scattering length by ex- 
ternal fields [[!^,[4(J. This allows a 'tuning' of the scattering length through both positive 
and negative values. Finally we recall that we have insisted on equal masses and trapping 
potentials for the two components. If these constraints are relaxed, we may more easily vary 
the values of the scaled interaction strengths. 

The oscillations in the fermion density distribution near the trap center reflect the matter 
wave modulation of the particles in the outermost shell. Their de Broglie wavelength can 
be estimated in the Thomas- Fermi approximation from (||): In the center of the trap the 
particles in the £ = states experience a vanishing potential for h = 0. As the Thomas- Fermi 
expression for the Fermi energy of Np fermions in a harmonic potential is Ep = (GNpY^hu 
we find for the de Broglie wavelength 

2tt \Z27ra 

DB ~ mo) ~ WW^ ~ a °' ( } 

an estimate that is reproduced by the data, see inset. 
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We now turn to the case of equal numbers of bosons and fermions. The influence of the 
inter-species interaction grows as the number of fermions is increased with dramatical effects 
on the atomic distibutions, as we shall demonstrate. We study the case of 10 6 fermions, and 
the same number of bosons with an interaction strength of g = 2.11/iu;ao 3 . We again expect 
that for certain critical parameters, the components find it energetically favorable to separate 
into two distinct phases, but this time bosons are expelled from the trap center, minimizing 
their internal interaction energy by spreading in a 'shell' around a fermionic bubble. Figs. |2| 
and |] present our results. The essential features are again the spatial separation of the two 
components, this time manifesting itself by the exclusion of the bosonic component from the 
trap center and the existence of a constant fermion density through the boson distribution 
for h = g. For a different choise of parameters, for example by letting the fermions be 
trapped by a weaker potential, we are also capable of producing a multi-layered structure 
with fermions residing on both sides of the bosons. 

We notice that as the bosons are expelled from the center of the trap, forming a 'mantle' 
around the fermions, the fermionic component is compressed, having a higher peak density 
and covering a smaller portion of the trapping volume. A similar behavior has been noted 
for bi-condensate systems |l£||Tj|. One of the essential features predicted in the Thomas- 



Fermi approximation is the existence of a 'plateau' of constant fermion density through the 
boson distribution for h = g. As illustrated by Fig. |3|, which is just a magnification of the 
central parts of Fig. |2|e, this phenomenon also appear in our quantum mechanical treatment, 
although with the parameters chosen it does not involve quite as many particles as obtained 



from the semi-classical calculations in 27 



It is interesting to compare the above mentioned results with those obtained by treating 
the fermions in the Thomas-Fermi approximation. This is done in Fig. £|for Nb = Nf = 10 6 , 
h = g = 2.11hual, and we note that the semi-classical description gives a qualitatively correct 
description, in that it reliably predicts the phase separation. Thus it is reasonable to use this 
approximative treatment of the fermions in the anisotropic case, where the exact description 
is too cumbersome. 



B. Anisotropic trap 

We now turn our attention to the anisotropic potentials, where we will use only the semi- 
classical Thomas-Fermi approximation for the fermion density. We aim to reveal similar 
variations in the ground state density profiles as for the isotropic trap, but going to higher 
dimensions we now have the opportunity to investigate the phenomenon of spatial symmetry 
breaking. Intuitively, we assume that for critical parameters it may be preferable for the 
two components to break mirror symmetry (z — > —z), thereby minimizing their mutual 
interaction, especially in elongated traps. Such a behavior has been predicted by Ohberg 
and Stenholm for bi-condensates in two dimensions ||18|| . 

It remains to be demonstrated though, that the features described in the case of the 
isotropic trap are still essential, when we consider the anisotropic scenario relevant in com- 
parison with currently experimentally feasible setups. We present in Figs. |5] and || the analog 
of Fig. |I] with the same choise of parameters and A = l/y/8, i.e. the inverse of the value 
for the current traps which have the strongest confinement along the z-axis. We notice 
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the appearance of the same qualitative features as in the isotropic trap, that is component 
separation for h > g and a plateau of constant fermion density for h = g. 

The 10 6 bosons are in the condensate which is unaffected in form and location by the 
presence of the relatively few fermions. Not shown in Figs. |5| and || is the distribution of 
fermions for h smaller than g. In this case the fermionic component overlaps the boson cloud 
at the center of the trap. 

To address the issue of symmetry breaking we adopt the same procedure as Ohberg and 



Stenholm ||18|| . This offers only suggestive evidence that symmetry breaking may occur. To 
investigate this behavior correctly one must use an altogether different approach, minimizing 
the energy functional to find the ground-state density profile P^| . The point is that the 
solutions of the Gross-Pitaevskii equation are stationary points of the energy functional, not 
necessarily corresponding to minima. They may therefore be unstable in certain parameter 
regions. It is possible though to single out the more stable of two configurations by comparing 
their total energy as this is minimum in equilibrium. 

The total energy functional of the two-component system is a sum of four terms 

E = T B + T F + V ext + V mt . (22) 

The first term is the boson kinetic energy 

Tb = Jd 3 r^\V^\ 2 . (23) 

As a fermion with wave number k(r) has a kinetic energy of h 2 k 2 /2M, the total fermionic 
contribution to the kinetic energy is found by integrating this local term over all of phase- 
space, weighted by the phase-space density, l/(27r) 3 , 

d 3 r r k F(r) h 2 k 2 
T P = / — — / d 3 k 



(2tt) 3 Jo 2M 

d 3 r h 2 r 9 ,^l5/3 , . 

— Q7r 2 n F (r) . 24 

2tt 2 10M L V ; J V ; 

Calculating the potential energy terms is easy, as they involve only integrals over the atomic 
densities 

V ext = J d 3 r X -Mu?r 2 [n B (r) + n F {r)} (25) 
V int = J d 3 rn B {f) \gn B {f) + hn F {f)} . (26) 

We have chosen the number of atoms to be iV# = N F = 10 3 , while the asymmetry parameter 
is still set to A = l/y/8. The interaction parameters are g = Q.Qlfiua^ and h = 5g. 

Starting the iteration with two well separated clouds displaced along the cylinder axis, i.e. 
along the direction of the weaker trapping potential, the calculation converges to a situation 
where the fermions localize on both sides of a central concentration of the Bose condensate: 
a 'boson-burger', see Fig. [7|. Initiating the calculation with two overlapping clouds in the 
center of the trap results in just the reversed situation: a 'fermion-burger', consisting of a 
central fermionic part surrounded on two sides by bosons, but this configuration has a larger 
energy. The 'boson-burger' seems to be the stable solution. 
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In Fig. [8] we show the spatial distribution of 5000 fermions and 1000 bosons. The par- 
ticles feel the same trapping potential as in Fig. |7], and the interaction strengths are kept 
unchanged. This configuration is the result when the starting point of the calculation is two 
separated clouds. When we start by placing both species at the center of the trap we achieve 
again the 'fermion-burger', but at a higher energy. Thus we conclude that in this region of 
parameter space the system is unstable against breaking of the reflection symmetry. 

We note that our approach provides two degenerate symmetry broken states, the one 
in Fig. [8] and its mirror image in the xy-plane. Going beyond our theoretical treatment 
(Hartree), we may construct superpositions of these two macroscopically states which do 
not break the spatial symmetry. One of these states will have a lower energy, but such a 



'Schrodinger-cat' state is exceedingly complicated to prepare, c.f. the discussion in [41 
Thus the symmetry broken solution is most likely to be observed in an experiment. 



V. CONCLUSION 

In this paper we have investigated the zero temperature ground state of a mixture of 
boson and fermion gases in both isotropic and anisotropic trapping potentials. We have 
addressed the issue of component separation using nummerical techniques to solve coupled 
equations for the spatial density of the two species. Our calculations have confirmed and 
expanded upon the results of a previous paper [^7| , which treated the problem only within the 



Thomas-Fermi approximation for both components and which analyzed only the case of an 
isotropic trap. We have confirmed the existence of three distinct states of the system under 
variation of the ratio of the interaction strengths h/g: For small values of this parameter the 
gases are interpenetrable, overlapping throughout the occupied volume of the trap, as their 
mutual repulsion is not strong enough to cause separation. When the coupling strength h 
exceeds the strength of the boson-boson interaction one of the species is expelled from the 
center of the trap. The spatial configuration in this case depends on the symmetry of the 
trapping potential. In an isotropic trap the separated phase is rotationally symmetric, the 
excluded component constitutes a spherical shell wrapped around a centrally compressed 
bulk. The anisotropic trap however has a parameter region where a breaking of symmetry 
(z — > —z) may occur, and we have demonstrated such forms. In the limiting case h = g 
there exists the possibility for the fermions to have a constant spatial density where the 
bosons are localized. 

An aspect of this work is the availability of an almost isolated degenerate Fermi gas 
through the complete separation of the two species. The trapped, degenerate Fermi gas is 
interesting in view of the possibility of a BCS transition when two spin states are trapped 
simultaneously |32| ,[42| , |43| and because of the analogies between this system and both atomic 



nuclei and the interior of neutron stars. 

The details of sympathetically cooling the Fermi gas to the degeneracy level through 
thermal contact with the Bose condensate are of course of great importance in further 



research ||44|| . In general the investigation of the cooling ability of the condensate should 
not be restricted to fermionic impurities. In view of the recent trapping of simple molecules 
in both optical ME\ and magnetic [ffij] potentials, also more complex solutes with several 



internal degrees of freedom pose an interesting challenge for future research. 
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Another direction worth noticing is the prospect of trapping a boson-fermion mixture in 
the periodic potential of an optical lattice [47], both in its own right and as a study of solid 
state phenomena. With quantum gases well beyond the degeneracy level a complete filling 
of the potential wells may well be expected |48[ . 

Finally it should be mentioned that in this work we have concentrated solely on systems 
with a positive coupling strength h. Allowing the interaction between the species to become 
attractive is known to induce a dramatic change in the macroscopic behavior of the system 
as it becomes unstable against collapse for large negative values of h [27j. We are cur- 
rently setting up calculations to investigate this phenomenon in detail using the numerical 
procedure developed in this work. 
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FIGURES 




FIG. 1. Spatial distribution of 1000 fermions as function of distance from the trap center. 
g = 0.015hu>aQ, and the boson- fermion coupling takes the values, h = (solid curve), h = g/2 
(long-dashed curve), h = g (dashed curve), h = 3g/2 (dotted curve), and h = 2g (dot-dashed 
curve). A magnification of part of the top curve is shown in the inset. The Bose-Einstein condensate 
component of 10 6 atoms is unaffected by the fermions and extends out to the distance ~ 7ao, where 
cusps in the fermion distributions are visible. 
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FIG. 2. Spatial distribution of 10 6 bosons (solid curves) and 10 6 fermions (dashed curves) as 
function of distance from the trap center. In Figs. 2a-2f, g = 2.llhtoaQ, and the boson- fermion 
coupling takes the values, h = 0, g/A, g/2, 3<?/4, g, hgj\. 
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FIG. 3. Spatial distribution of 10 6 bosons (solid curves) and 10 6 fermions (dashed curves) as 
function of distance to the trap center, h = g = 2.1\Tiuj(1q 
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FIG. 4. Comparison of density profiles calculated by using the Slater-determinant method 
(solid lines) and the Thomas-Fermi approximation (dashed lines) for the fermion density. 
N B = N F = 10 6 , h = g = 2.11/iwag. 
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FIG. 5. Spatial distribution of 1000 fermions in an anisotropic trap with A = 
g = 0.015/iu;ao, and the boson-fermion coupling takes the values, h = g. The Bose-Einstein 
condensate component of 10 6 atoms is unaffected by the fermions and extends out to the distance 
r ~ 7ao- 
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FIG. 6. Same as Fig. ||, but with h = 2g. 
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FIG. 7. The 'boson-burger' configuration of 1000 bosons (solid lines) and 1000 fermions (dashed 
lines) in a prolate harmonic trap with A = 1/ y/8. The strength of the interparticle interactions are 
g = 6.67htoaQ and h = 5g. The iteration was started with two Gaussian density profiles located 
oppposite to each other, away from the center. 
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FIG. 8. Spatial density of 5000 fermions and 1000 bosons in a prolate trap with A = l/\/8 and 
interaction strengths g = G.67hu}a,Q and h = 5g. The bosons (solid lines) and the fermions (dashed 
lines) localize in separate regions of space due to their strong repulsion and the elongated nature 
of the trap. 
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